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In this work we perform the parameter tuning of three flavors of dynamical clover quarks on 
anisotropic lattices. The fermion action uses three-dimensional spatial stout-link smearing. The 
gauge anisotropy is determined in a small box with Schrodinger background using Wilson-loop ratios. 
The fermion anisotropy is obtained from studying the meson dispersion relation with antiperiodic 
boundary conditions in the time direction. The spatial and temporal clover coefficients are fixed to 
the tree-level tadpole-improved values, and we demonstrate that they satisfy the nonperturbative 
conditions as determined by the Schrodinger-functional method. For the desired lattice spacing 
a s ~ 0.12 fm and renormalized anisotropy £ = 3.5, we find the gauge and fermionic anisotropies can 
be fixed to quark mass independent values up through the strange quark mass. This work lays the 
foundation needed for further studies of the excited-state hadron spectrum. 

PACS numbers: 11.15.Ha,12.38.Gc,12.38.Lg 



!>• ' I. INTRODUCTION 

(N ' 

Lattice quantum chromodynamics (QCD) has successfully calculated many properties of the hadronic spectrum. 
However, there remain many challenges for the lattice community to resolve in determining the myriad states present 
i—* \ in QCD. The case of the nucieon spectrum is one such example. Consider the lowest three states in the TV spectrum 
O ■ (TV, TV' (Pn) and TV* (Sii)), for example. Many earlier quenched lattice QCD calculations [3, H, H, 0, H, @1 find a 
O | spectrum inverted with respect to experiment, with the TV' heavier than the opposite-parity state TV*. Although the 
i-C ■ Kentucky group Q managed to find the correct mass ordering around a pion mass of 300-400 MeV (after taking care 
of the effects of the quenched "ghosts" ) , no other lattice group has been able to reproduce the experimental ordering 
, using different approaches. Furthermore, these are just the lowest few states in the TV spectrum. There are many 
' more states seen in experiment for which lattice calculations could help in the identification of particle properties. 

This situation suggests an urgent need for full-QCD simulations that can resolve some of these issues. In order 
to improve the signal of the excited states (especially for the higher-excited nucieon spectrum), one needs a lattice 
with a fine temporal inverse lattice spacing on the order of 6 GeV. At the same time, we also want to avoid finite- 
volume effects. Current dynamical lattice gauge ensembles manage to have a reasonable lattice volume with spatial 
dimensions of about 3 fm; however, typically the inverse lattice spacing is about 2 GeV, which is not fine enough 
to allow an accurate determination of more than one excited state. One solution to this problem is to generate 
anisotropic dynamical lattices. 

Anisotropic techniques have been widely adopted in lattice calculations. Anisotropic rclativistic heavy quark actions 
have been used for charmonium studi es la, B . Another main application is for calculations, such as glueballs [lfjl ] and 
multiple excited-state extraction [TTL Il2l. Il3l| , where the anisotropic lattice technique has advantages over isotropic 
lattices due to the finer temporal lattice spacing. However, there is a worry that uncontrollable quark mass effects, 
ma s , might enter into systematic errors [LH flH Il6l] when m (the quark mass) is large and the spatial lattice spacing 
a s is about 0.1 fm. Since we are working in the light-quark limit, this is not a major concern here. In order to 
remove possible 0(a) systematic errors from the action, we tune our anisotropic action to be as close to the (on-shell) 
nonperturbatively correct action as possible. 

Previous results on anisotropic lattices include two-flavor anisotropic dynamical simulations done by CP-PACS [ij} 
and TrinLat collaboration CP-PACS performed the first dynamical TV/ = 2 anisotropic tuning of dynamical 

clover fermions [T3] (without gauge- link smearing), using the Iwasaki gauge and Sheikholeslami-Wohlert (also called 
clover [lj|) fermion actions. In that study, they set the coefficient of the clover term within the clover action to 
tadpole-improved tree-level values. The TrinLat collaboration [l8j used a two-plaquette Symanzik- improved gauge 
action with tree-level tadpole improvement and a Wilson fermion action with a Hambcr-Wu term. One should also 



'Electronic address: cdwards@jlab.org 
T Electronic address: bjoo@jlab.org 
i Electronic address: ihwlin@jlab.org 



2 



note that they adopted stout-link smearing [20| of the spatial gauge fields in the clover action. That is. the gauge 
fields entering the fermion action were not smeared in the time direction, preserving the positivity of the fermion 
transfer matrix. Only two iterations of stout smearing were used with a staple weight p = 0.22. 

In this work, we will use a three- flavor clover action with stout-link smearing (in the spatial directions only), and an 
0(a 2 )-improved Symanzik gauge action. Working in the Schrodinger- functional scheme [U HU H3 , we determine 
the gauge anisotropy by computing Wilson loop ratios with the background field applied in the z direction. The 
fermion anisotropy is determined from the conventional meson dispersion relation with periodic boundaries in the 
spatial directions and antipcriodic boundaries in the time direction. The coefficients of the gauge action are set to the 
tree-level tadpole-improved values, and the clover coefficients are fixed at the tree-level stout-link smeared tadpole- 
improved values, with the tadpole factors set from numerical simulation. We demonstrate that the clover coefficients 
are consistent with nonperturbative values determined in the Schrodinger-functional scheme. Our configurations have 
been generated using the Chroma [2{| HMC code with RHMC for all three flavors and multi-timescale integration. A 
preliminary study can be found in Ref. [26| . 

The structure of this paper is as follows: In Sec. [ill we will discuss the details of the actions used in this work, the 
stout-link smearing and Schrodinger-functional scheme calculations, and how we determine the coefficients. Then we 
will cover the Rational Hybrid Monte Carlo (RHMC) used in this work, how we apply it on anisotropic lattices with 
even-odd preconditioning, and how to use these techniques with stout-link smearing in Sec. lIIIl We present numerical 
results in Sec. IIV1 where the gauge and fermion anisotropy, and PCAC mass are measured, and their corresponding 
(tuned) bare values are determined. Some conclusions and future outlook are presented in Sec. [V] 



II. METHODOLOGY AND SETUP 



A. Action 



In this section, we describe the gauge and fermion actions used in this calculation. For the gauge sector, we use 
a Symanzik- improved action which was used in the glucball study of Ref. With tree- level tadpole-improved 

coefficients, the action is 
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where V is the plaquette and 1Z is the 2x1 rectangular Wilson loop. The coupling g 2 appears in (5 = 2N c /g 2 . The 
parameter £ is the bare gauge anisotropy, and u s and u t are the spatial and temporal tadpole factors, dividing the 
spatial and temporal gauge links, respectively. This action has leading discretizations error of 0(«s, a 2 , g 2 a 2 ) and 
possesses a positive definite transfer matrix since there is no length-two rectangle in time. 
In the fermion sector, we adopt the anisotropic clover fermion action Q 
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In terms of dimensionless variables ip = a^ 2 ip, rho = nioa tl V M = a^V^, A M = a^A^, = a^a^F^ and the 
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dimensionlcss "Wilson operator" = V M — ^7 M A M , we find the fermion matrix Q becomes 
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Here v is the ratio of the bare fermion to the bare gauge anisotropy. From the field redefinition [27], |28| , there is one 
redundant coefficient: cither v t or v s . There arc two common choices to eliminate this redundancy: setting v s = I 
(i/(-tuning) or v t — 1 (iAj-tuning) [sL We will use ^-tuning in this work, so we set vt = 1, with the tree-level conditions 
on c s and c t as described in Ref. [S]. In particular, we choose the tree-level tadpole- improved values 
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where u s and itt are the tadpole factors and the fraction at/a s = 1/C is set to the desired renormalized gauge 
anisotropy. 

In this work, the gauge links in the fermion action are 3-dimensionally stout-link smeared gauge fields with smearing 
weight p and n p iterations. To distinguish tadpole factors associated with the smeared fields appearing in the fermion 
action from those appearing in the gauge action, we use notations u s and Ut for spatial and temporal tadpole factors, 
respectively. For convenience of parameterization, we use the bare gauge and fermion anisotropies, J g j, defined as 



7 9 = Co, 

To summarize, the final gauge and fermion actions are 
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B. Stout-smeared links 



The smeared fermion action provides significant improvements on actions that explicitly break chirality such as the 
clover fermion action. It has been demonstrated that chiral symmetry is improved [29j | after treatment of the gauge 
links in the fermion action. In this work, we use three-dimensionally stout-smeared links [20| in the fermion action. 
We use smearing parameters p = 0.22 and n p = 2 (as in Ref. [lH) through the end of SecHU In the numerical section, 
we will examine our choice of stout-smearing parameters with greater caution. Since the smearing does not involve 
the time direction, the transfer matrix remains physical. As with other smearing techniques, we need to check the 
smearing parameters carefully to avoid potentially incorrect short-distance physics. 

We consider the effects of our choice of action on scaling 30] in a quenched theory. On the left-hand side of Figure [fl 
we show the scaling behavior of the vector meson mass in units of the string tension, amy /\/a 2 a. In this plot, all 
the points on the graph have fixed quark mass determined from the ratio of the pseudoscalar to vector meson mass 
ratio mps/rnv = 0.7. The left-hand panel shows results for the quenched Wilson gauge and fermion action, and 
the right panel shows results for the quenched anisotropic Wilson gauge and clover fermion actions. The curves are 
scaling fits to the Wilson and clover fermion data constrained to have the same continuum limit [3(| , with the Wilson 
action scaling like 0(a) and the non-perturbative clover results scaling like 0(a 2 ). The horizontal line is the (fitted) 
continuum limit value. The small residual scaling violations in the nonperturbatively improved clover action (1% at 
a 2 <j ~ 0.05 or a ~ 0.1 fm) indicate that the dominant source of scaling violations in the Wilson action comes from 
chiral symmetry breaking at 0(a). Our simulations on anisotropic lattices show similar scaling for the Wilson fermion 
action (the diamonds), and when we add stout-link smearing in the fermions (the squares), we see a large reduction in 
the scaling violations. In the clover case (on the right-hand side of Figure[T]), similar tests were performed, with both 
tadpole-improved perturbatively determined and nonperturbatively determined clover coefficients. When stout-link 
smearing is added to the anisotropic clover action, the scaling remains consistent with the unsmeared results. Thus, 
the three-dimensional stout-link smearing does not adversely effect the scaling properties of the fermion action, and 
the resulting scaling violations after smearing arc suitably small. 
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FIG. 1: The scaling behavior of the quenched Wilson gauge action with Wilson (left) and clover (right) fermion actions showing 
the effects of anisotropy and stout-link smearing. Shown are scaling fits of the Wilson and clover results constrained to have 
the same continuum limit [3(j. In the left plot, the anisotropic stout-link smeared result is the square. In the right plot, the 
"NP" labels non-perturbatively tuned clover coefficients, and the "Tad" indicates tadpole-improved coefficients. 



C. Schrodinger functional 

The Schrodinger-functional scheme [111, [12, allows for simulations at small pion mass, since the background 
field lifts zero modes. We use the Schrodinger functional and the PCAC relation to check how close our csw m the 
fermion action is to the nonpcrturbative value. In previous work with dynamical fermions, the Alpha collaboration 
used two-flavors [HI with the Wilson gauge a ction, and CP-PACS used two-flavor and three-flavor calculations with 
the Wilson and Iwasaki gauge actions [171 |32|. |33| . In this work, for our calculations with the background field in the 
spatial direction where we have length-two rectangles in the gauge action, we used the Schrodinger-functional setup 
from Ref. [H]. 

The bare PCAC quark mass is calculated through the PCAC relation using zero momentum projected correlators 

as 

m(x ) = r(x ) + ac A s(x ), (9) 

where xq is some time slice (or possibly space slice) of the correlator, ca is the O(a)-improvod coefficient for the axial 
current and 

1 (S + d$) /aOo) 
4 fp{xo) 

1 d d^fp(x ) , . 
s(x ) = -a — - - - . (10) 

2 fp[xo) 

The axial current j 'a with F = 7.57^ (or pseudoscalar density fp with T = 75) is a correlation function of bulk fields 
(-0, "0) and boundary fields at t = (77, rf) defined as 

x y,z 

Similarly, correlators propagating from the other wall lead to definitions of r' and s', where the f A r P \ now involves 
the other boundary fields at t = T (jf, 77') and a sign change. 

On an isotropic lattice, one can also determine the clover coefficient csw in Schrodinger-functional scheme from the 
PCAC relation. For an arbitrary set of action parameters, the relation m(yo) — m'(yo) is not satisfied in general is 
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not satisfied for a generic j/o- We define an intermediate (before nonperturbativc tuning) ca as ca according to 



- / x 1 r'(y ) - r(y ) 
a s'{y ) - s{y ) 



(12) 



and define a modified mass in Eq. [9] as 



M(x ,y ) = r(x ) - c A (yo)s{x ) 
M'(x ,yo) = r'(x Q ) - c A (y )s'(x ). 

One nonperturbatively determines cqw by imposing the condition 

AM = M(x , yo) - M'(x ,y ) = AM' ', 



(13) 



(14) 



where AM' ' is either set to be zero or the tree-level mass splitting which can be obtained from a free-field simulation 
with the same setup of the gauge and fermion actions. Note that the choice of (xo,yo) is (T/2, T/4) for both AM 
and M. csw is obtained by tuning the condition AM' ' — AM = 0, and ca is obtained in straightforward fashion 
from Eq. [T2] with the correct csw- Note that all applications of the Schrodinger functional so far have been limited 
to isotropic lattices. 
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FIG. 2: Amin(Q 2 ) measured in three-flavor anisotropic clover simulations (using parameters (3 — 2.2, p — 0.22, n p = 2, 
mo = —0.054673, 7/ = £ = 3.5) with (left) and without (right) a background field. The i-axis is in units of 5 trajectories. The 
right-hand panel (no background field) has a longer autocorrelation time compared to the case with a background field (left 
side). The vertical scale, however, is different because of the background field. 
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FIG. 3: a s Ms comparison among different spatial volumes. The parameters in the action are the same as Figure [2] The right 
panel is an enlarged scale of the left panel with the addition of the 12 3 spatial-volume result, which is the size used in the 
remainder of this work. 
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In this work, wo implement the Schrodinger-functional setup on anisotropic lattices for the first time in dynamical 
simulations. In one of our earlier three- flavor (anisotropic) clover action simulations (with parameters f3 = 2.2, 
p = 0.22, rip = 2, mo = —0.054673, 7/ = £ = 3.5), we found that the autocorrelation time of the lowest eigenvalue of 
Q^Q (as defined in Eq.[3J is significantly reduced by introduction of the background field, as shown in Figure[5J Thus, 
for these background-field calculations, we can use correspondingly fewer trajectories compared to the calculations 
with antiperiodic boundary conditions. 

We also implement the background field not only in the t direction (as conventionally used in Schrodingcr functional) 
but also the z direction. (We label the modified masses M s and M t with respect to the direction of the boundary 
field direction and similarly for the mass difference AM S and AM t .) We measure M(xo) with various spatial x and 
y sizes, L x y , and sizes L Z1 showing results in Figure [31 all with time extent L t = 32 and equal amount of statistics. 
When we increase the length in the z direction, a good signal appears from t = 12 to 16 but not beyond 24. This is 
because the background field becomes too weak at large L z . When we increase the background field signal, which is 
proportional to $ at L z = 24, the signal is still not as good as for L z — 16. Similar checks can be done regarding the 
size of L x . y . As we increase the value from 8 to 12, the signal shows improvement. The right panel of Figure [3] zooms 
in the scale as seen from the left panel, and demonstrates that the L x . y = L z = 12 spatial volume has the lowest 
statistical error while also providing enough spatial size to adequately resolve the potential which will be discussed in 
Section [IV Bl this is the spatial volume we use in the remainder of this work. 

We note that the fermion algorithm we will use corresponds to simulating the Nf = 3 fermion determinant as 

3 

det (Q^Q) 2 ■ In particular, when the renormalized quark mass is negative (below the chiral limit), there is no phase 
that should arise as in a Nj = 3 version of QCD. Thus, the algorithm is implementing the absolute value of the 
fermion determinant. While, strictly speaking, this is not a Nf = 3-flavor version of QCD when the (renormalized) 
quark mass is negative, it does allow us to implement a SU(3) flavor symmetric version of PCAC. 



D. Renormalization conditions 



In principle, we can determine the critical values for the bare parameters 7* and all the unknown coefficients in the 
fermion action nonperturbatively by imposing the following conditions: 

€g(lg,lf,™0, = i 

Mt(l* g ,l},rn* ,c* s ,c* t ) = m q 
AM s (7 s *,7/,m5,c:,c t *) = AM<°> 

AM t (7 s *,7/,mS,c:,c t *) = AM t (0) , (15) 

where £ ff is the renormalized gauge anisotropy, £/ is the renormalized fermion anisotropy (defined through the meson 
dispersion relation E(p) 2 = M 2 +j? 2 /£/) an d M t is the PCAC quark mass measured in Schrodinger-functional 
background field (from Eq. [5]). AA/ S t are the mass differences measured in two different background field directions. 
However, this requires that we search a five-dimensional parameter-space with dynamical simulations which is costly. 

In this work, we set the spatial and temporal clover coefficients c s and Ct (in Eq. [5]) to the tree- level tadpole-improved 
values. (Later in the numerical section, we demonstrate using these values that the last two conditions in Eq. 1151 do 
hold within a few percent.) To determine the remaining three coefficients, we parameterize £ s , and M t as functions 
fiilgi 1ft m o)- For simplicity, we choose functions of the bare parameters with the form 

£g(lg>lf, m o) = ao + ailg + a2lf + a 3 m 
€f(lg,lf> m o) = b o + hlg + hlf + b 3 m 

M t (ig 1 7/i m o) = c + ci7 9 + c 2 7/ + c 3 m . (16) 

These parameterizations are linear functions in the coefficients. We can choose higher powers of the bare parameters 
such as a TOq term; however, the coefficients that we need to determine would remain linear. 
Once the coefficients Oj, 6^, and are determined, we impose our renormalization conditions 

€g{lg,lf,mo) = £ 
Hf(l* g ,lf,rn* ) = € 

M t (l* g ,lf,m* ) = m q (17) 
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to obtain the critical values for the bare parameters as a function of the input quark mass ra q . If only linear terms in 
the bare parameters are used, then the intersection of these three hyperplanes is the solution of a 3 x 3 linear system 
of equations. If higher-order terms are used, then the intersection is the root of a system of functions. The fitted 
parameters determine the chiral limit when the input m q = 0. 



III. ALGORITHM 

Our configurations were generated with the Rational Hybrid Monte Carlo (RHMC) algorithm [H, HE Hfil • Strictly, 
RHMC refers only to the method for simulating odd flavors of fermions, and one can combine several orthogonal 
algorithmic improvements with the RHMC scheme resulting in a wide variety of possible RHMC algorithms. In order 
to be specific therefore, we describe the RHMC method in brief below and then detail our particular combination of 
improvements. 



A. Rational Hybrid Monte Carlo 

The basic technique for gauge generation is a Markov Chain Monte Carlo method, where one moves from an 
initial gauge configuration to a successive one by generating a new trial configuration and then performing an ac- 
ceptance/rejection test upon it. If the trial configuration is accepted, it becomes the successive configuration in the 
chain, otherwise, the original configuration becomes the next state in the chain. 

In order to use a global Metropolis accept/reject step with a reasonable acceptance rate, the space of states is 
extended to include momenta 7r^(x) canonical to the gauge links U^(x) so that one may define a Hamiltonian 

h = \Y,^ x )^^ x ) + s ( u ) ( 18 ) 

where S is the action. It is then possible to propose new configurations from previous ones by performing Molecular 
Dynamics (MD). Using a reversible and area preserving MD step maintains detailed balance, which is sufficient for 
the algorithm to converge. In order to ensure ergodicity in the entire phase space, the momenta need to change 
periodically. This can be affected by refreshing the momenta from a Gaussian heatbath prior to the MD update step. 

In order to deal with the fermion determinant, it is standard to use the method of pseudofermions. One integrates 
out the Grassman- valued fermion fields in the action and rewrites the resulting determinant as an integral over bosonic 
fields 

Z = [ [dfj}[dr]}e- rrD ^ = det (X>) = / [d^p^e-^^ (19) 



where r\ and fj are the Grassman valued fields, D is some fermionic kernel and <jr and <f) are the bosonic pseudofermion 
fields. Our phase space is thus enlarged to include also the pseudofermion fields, which similarly to the momenta, 
need to be refreshed before each MD step, to carry out the the pseudofermion integral. 
In the case of a two-flavor simulation, T> is typically of the form 

V = Q f Q (20) 

where Q is the fermion matrix for an individual flavor of fermion. In this case T> is manifestly Hermitian and positive 
definite, and the integral in Eq.[l9]is guaranteed to exist. Furthermore, the pseudofermion fields can easily be refreshed 
by producing a vector x filled with Gaussian noise with a variance of | and then forming = Q'X- 

In th e case of an odd number of flavors, since Q itself is not guaranteed to be positive definite, one works instead 
with y/WQ)- 

det (Q) = det (Q f Q)^ - f [d^H^e"^^"^ . (21) 



The square root in T> can be approximated to numerical precision using a low-order rational approximation r 5 (Q^Q) 
which can be expressed in a partial fraction (sum over poles) form as 

(QtQ)-^ r -|(QtQ ) = ao/ + ^ a ^QtQ + ^]- 1 , (22) 
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where the rational function approximation is specified by the coefficients ctk and /3k- In particular the action of all 
the denominator pieces onto one vector <p involves the solution of linear systems 

[Q^Q + (3k}x = <t> (23) 

which can be performed simultaneously using a multiple-shift (a.k.a. multi-mass) conjugate-gradient solver [37| . 
Refreshment of the </> field now proceeds by evaluating 

<j> = A(Q^Q)r 1 (24) 

. i 
where the rj are once again filled with Gaussian noise of variance \ , and ri is now a rational approximation to (Q'Q) 4 . 
Molecular dynamics forces now need to be calculated for each pole term in the partial fraction 

National = - ^ °* ^ (^Q + #0 
fc 

The idea of using a rational approximation in the action, and its consequences for field refreshment and molecular 
dynamics forces make up the basics of the RHMC algorithm. 

In terms of tuning, in our single-precision simulations we tune the rational approximation coefficients by requiring 
that the approximation has a maximum error over the approximation interval that is smaller than the solver residua 
we require in the evaluation of the partial fraction expansions, namely 10 -8 for energy calculations and 10~ 6 for the 
force calculations. Further, in the case of the force calculations we successively relax the solver criteria for poles that 
have smaller contributions to the force as in Refs. [34l. l35l. |3(| . 



(Q f Q + fa) <t>. 



(25) 



B. Multi-Scale Anisotropic Molecular Dynamics Update 

While any reversible and area-preserving MD update scheme can be used in the MD step, the acceptance rate is 
controlled by the truncation error in the scheme. This manifests itself as a change in the Hamiltonian, SH, over an 
MD trajectory, since we use the Metropolis acceptance probability 

P acc =min(l,e- 5ff ) . (26) 

We may easily construct a manifestly reversible scheme by combining symplcctic update steps U p (St) and U q (5r) 
which update momenta and coordinates by a time step of length St respectively 

Uri&rJ: {tt^U^x)) -» (tt m (i) + F M (a:)*r MJ U^x)) (27) 
U q (5T^): (tt^U^x)) -> (ir^x),e Mr »U^x)) , (28) 

where F^i^x) is the MD force coming from the variation of the action with respect to the gauge fields. We emphasize 
that one may update all the links pointing in direction /i with a separate step-size St^. While this may not be useful 
in isotropic simulations, in an anisotropic calculation with one fine direction, it may be advantageous to use a shorter 
timestep to update the links in that direction to ameliorate the typically larger forces that result from the shorter 
lattice spacing [l8| . The anisotropy in step size requires a small amount of manual fine tuning, but should be similar 
to the anisotropy in the lattice spacings. 

Our base integration scheme in this work is due to Omelyan (38l . l39l . |40| ; we use the combined update operator 

U^St) = Kp(XST)l( q (-5T)Up(l - 2X5T)U q {^5T)U p {XST) (29) 

which results in a manifestly reversible scheme that is accurate to 0(<5r 3 ). The size of the leading error term can be 
further minimized by tuning the parameter A. In our work we used the value of A from Ref. [39l | without any further 
tuning, which promises an efficiency increase of approximately 50% over the simple leapfrog algorithm. 

In Refs. fiol . |4]| it was shown that a reversible, multi-level integration scheme can be constructed which allows 
various pieces of the Hamiltonian to be integrated at different timcscalcs. Let us consider a Hamiltonian of the form 



H(tt, U) = -Trj, (x)7r„(a;) + S 1 {U) + S 2 (U) 



(30) 
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where S\(U) and 52(f) are pieces of the action with corresponding MD forces F\ and F 2 respectively. One can then 
split the integration into 2 timescales. One can integrate with respect to action Si(U) using U 1 (5ti), where in the 
component U p (Sti) we use only the force Fx. The whole system can then be integrated with the update 

U 2 (St 2 ) =U' p {\8t 2 )U 1 {^5t 2 )U' p {\ - 2\5t 2 )U 1 {^8t 2 )U' p {\8t 2 ) (31) 

where in U' p we update the momenta using only F 2 . Thus we end up with two characteristic integration timescales 
5tx and St 2 . The scheme generalizes recursively to a larger number of scales. A criterion for tuning the algorithm 
is to arrange for terms in the action to be mapped to different timescales so that on two timescales i and j we have 
I (Filler, « H-Fjll^Tj, as suggested in Ref. [Hf. We now proceed to outline how we split our action. 



C. Even-Odd Preconditioning the Clover Term 

The clover term may be preconditioned by labeling sites as even and odd, and grouping together the terms in the 
operator that connect sites of various labels. In particular 

« - ( t aZ ) ■ (32 » 

where the blocks A ee and A ao contain the clover term and the diagonal parts of the Wilson operator whereas the D eo 
and D oe contain the Wilson hopping term. This matrix can be block diagonalized as 



1 \ / A ee W 1 A-JD 



D oe A~} 1)1 A ao - D oe A~}D eo J \ 



(33) 



and it is clear that 

det (Q) = det (A ee ) det (A oa - D oe A^D eo ) . (34) 
We can write an action containing Nf degenerate fermion flavors as 

det(Qt Q )^ =e ^^ 1 o g A„(,)- 2: ri 1 ,t[ r -i (( 5t ( 5 ) ] 0i (35) 

with the preconditioned fcrmion matrix 

Q = A ao - D oe A~*D eo . (36) 

This manner of preconditioning the clover term is quite standard; guides to implementation are detailed in Ref. [42l |. 
We mention that it may be possible to combine the three pseudofermion terms each containing r~ 2 (Q^ Q) into a single 
one containing instead a single rational approximation r~5 ; however the Remez algorithm for this approximations 
results in negative roots with dire implications for our multi-shift solver. In this work therefore, we have simulated 
with three separate one-flavor pseudofermion terms. 

It is our experience that the forces resulting from the tr log A ee terms were at least an order of magnitude smaller 
than the pseudofermion terms, and so the two kind of terms could be run on separate timescales. 



D. Stout-Link Smearing in Fermion Actions 

The fermionic terms in our action employ stout-link smearing [20| on the links in the spatial direction. We leave the 
temporal direction unsmeared to keep the transfer matrix physical. The Schrodingcr-functional boundary condition 
is imposed at every iteration of the stout-link smearing. Our fermion operator Q is evaluated on the stout-smeared 
fields. In the calculation of the fermion forces, we compute the force on the stouted links, but then have to also apply 
the chain rule to compute the force coming from the original thin links 

dQ _ dQ dU stout 

dUthin dt/stout dUthin 

In particular the t ^ t °" t term is common to all the poles in the rational function approximation; hence, in our rational 
forces computations we compute all the forces with respect to the stout links first and then perform the recursion to 
the thin links only once. The Schrodingcr-functional boundary conditions imply that the force for the links which we 
hold fixed are set to zero. 
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E. Split Gauge Term 

We can write our gauge action schematically as 

S = S s (U) + S t (U), (38) 

where the term S s contains only loops with spatial gauge links, and the St term contains loops contain spatial and 
temporal links. While the term S s produces forces only in the spatial directions, the St term produces forces in both 
the spatial and the temporal directions. In particular the spatial forces from St are larger in magnitude than the 
spatial forces from S s by roughly the order of the anisotropy, and in turn, the temporal forces from St are larger 
than the spatial forces from St. Our anisotropic integration step-size balances the spatial and temporal forces of the 
St term against each other. However, in order to balance the spatial forces from St and S s against each other, we 
integrate them on separate time scales. 

F. Summary 

To summarize, we use the RHMC algorithm, with approximations tuned separately for the force and energy calcu- 
lations. Our molecular dynamics scheme uses anisotropic timcsteps and a recursively defined multi-level integration 
scheme based on Omelyan's inverter. Our fermionic terms employ stout smearing in the spatial directions only. We 
have 4 kinds of terms in our molecular dynamics integration: the tr log A ee term, the 0jr~2(Q)(/> pseudofermionic 
terms, the spatial gauge action term S s and the temporal gauge action term St- 

Based upon the magnitude of the molecular dynamics forces, we split our integration scheme onto three time scales 

• Time scale 1 uses the Omelyan integrator with a timestep of 6t±, and contains the tr logA ee term and the 
pseudofermion terms. 

• Time scale 2 uses a leapfrog integrator with a timestep of Sti relative to time scale 1, and contains the spatial 
gauge action term. 

• Time scale 3 uses a leapfrog integrator with a timestep of St 3 relative to time scale 2, and contains the temporal 
gauge action term. 

For our 12 3 x 32 volume results with a background field, we used timesteps {Sti, 5t2, St^) = (1/4, 1/4, 1/3) for the three 
time scales, and for our 12 3 x 96 volume results with antiperiodic boundary conditions we used (1/5,1/3,1/2). In 
addition, the time step for the temporal direction was a factor of £ = 3.5 times smaller than the spatial timcsteps. The 
acceptance rate was typically between 60 and 70%. These technologies arc all implemented in the Chroma software 
system (25| . 

IV. NUMERICAL RESULTS 

We are interested in a spatial lattice spacing on the order of 0.1-0.2 fm with a target anisotropy of £ = 3.5, which 
would provide a fine-enough temporal lattice spacing for excited-state physics. We proceeded by making an initial 
guess for the anisotropy parameters, "f g = 7/ = £ = 3.5 and mo = —0.05, and computed the Sommer scale ra/a s 
(Ref. [HI) over a range of (3 with lattice sizes of 12 3 x 32. Anticipating significant running in ro/a s , we chose to use 
f3 = 1.5. A preliminary investigation indicates a lattice spacing of roughly 0.12 fm. However, a careful determination 
of the lattice spacing involves a study of the static quark potential as well as hadron masses on Nf = 2 + 1 ensembles, 
a determination of the strange quark mass, extrapolation of the light quark mass to the physical limit, and examining 
different approaches to include potential systematic uncertainties which will be presented in Ref. [44| ■ Note that in the 
earlier part of this section, we will vary (3 to make a global search and to get an understanding of the stout-smearing 
parameters and tadpole factors. 

We tune the gauge and fcrmion anisotropics in three- flavor simulations. The determination of the tadpole factors 
used in the gauge and fermion action and the stout-link parameter study will be described in Sec. IIV Al Wc ultimately 
decide to 3d stout-link smear with n p = 2 and p = 0.14. Sec. IIV Bl and IIV CI describe how the gauge and fermion 
anisotropies are measured, and finally the target coefficients in the fermion action are determined in Sec. IIV Dl 
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A. Plaquette, Tadpole Factors and Stout Links 



In this work, we tadpole improve [4a | the gauge and fermion actions. This procedure amounts to replacing the 
gauge link fields with U p (x) = U p (x)/u. We obtain the tadpole factors used in the simulation from a nonperturbative 
tuning. There are two types of tadpole factors; we denote using u Sj t the tadpole factors with the stout-smeared fields 
(which are used in the gauge links of the fermion action), and u Stt are those without smearing (which are in the gauge 
action). For each parameter set, we start with a pcrturbative guess, u s = u t = 1, and calculate the actual u s-t from 
the plaquette. Then, we use the new value to feed back into the next dynamical run. The value of u s j will converge, 
giving the nonperturbative value with 1-2% precision. 

We calculated these tadpole factors with a large range of (3 values in three-flavor simulations. Over a wide range of 
parameters, we find that the square root of the spatial plaquette is consistent with the temporal plaquette, for both 
smeared and unsmeared links, indicating that the temporal tadpole factor is close to one. We chose to fix it 4 = u t = 1 
throughout the remainder of the work. It still remains to determine the spatial tadpole factors. 

Let us first concentrate on the influence of the stout-smearing parameters on the plaquette. In Sec. [Til we se t the 
stout-link parameters to be p = 0.22 and n p = 2, following the choice of the two-flavor anisotropic study done in 
Rcf. [18j . At one-loop level we expect the spatial plaquette to be 

(V s ) = (39) 

where cis' 1 is a function of the gauge action and stout-link parameters. Figure [4] shows ci)? as a function of p for 
rip = 1, 2, 3, as calculated in Ref. [46| with our choice of a Symanzik-improvcd gauge action. We find that the choice 

of p = 0.14 corresponds to a minimum of c S s which corresponds to the maximum of the spatial plaquette in Eq. [39] 
for g 2 > 0. We note that this value of p = 0.14 is also consistent with the suggestion [47| that the maximum allowable 
value, given by a classical argument, is 1/6 in our case of three-dimensional smearing. We also investigate the spatial 
plaquettcs in our numerical studies; see Figure [S] In this investigation, we use three different three-flavor ensembles 
with sea-sector parameters as follow: 



Ensemble A: (3 = 2.0, -y g = 3.5, 7/ = 3.89, {p,n p } = {0.22,2} 

Ensemble B: (3 = 2.0, j g = 3.5, 7/ = 3.5, {p, n p } = {0.14, 1} 

• Ensemble C: f3 = 1.0, j g = 3.0, 7/ = 3.33, {p,n p } = {0.22,2} 

with fixed niQ = 0, 12 3 x 32 volume with antiperiodic boundary conditions in time. We vary p in the valence sector 
and study the resulting behavior of the spatial plaquettes. We find that plaquette is maximized in the vicinity close 
to our initial choice p = 0.22. 

We then concentrate on one set of fermion coefficients, p = 0.14, uiq = 0, 7 S = 3.5, and 7/ = 3.5, and study the 
tadpole factors as functions of (3 with two values of n p , as shown in Figure [U We see that the fermion stout-link 

smeared tadpole factor u s (denoted as u s stS> in the plot) increases when we double the number of smearing steps and 
gets closer to one, which is what we expected. The gauge tadpole factor (without stout-link smearing, denoted as 

u£ ) remains unchanged, as if the stout-smearing in fermion sector had no impact in the gauge sector. This tells us 
the ultraviolet observable does not change with n p = 2; therefore, in the rest of this paper, we will fix the smearing 
parameters to p = 0.14 and n p = 2. 

We parameterize the smeared and unsmeared spatial tadpole factors as 



3 1 , 

E l + a n i 



2n 

(40) 



,2/i 



with the constraint a\ — b\ = c^, where g 2 is defined as 6//3 and cS 1 ^ is the one-loop perturbative value of the 
tadpole factor [46]. Figure [7] shows the data for n p = 2 and p = 0.14 and the fit using Eq. 201 The tadpole factors 
are interpolated well over a large range of /3. At j3 = 1.5 (on which we will focus for the rest of the paper) this 
parameterization gives: 

u s = 0.7336, u s = 0.9267. (41) 



In all further investigations, we use the values as determined by Eq. 0D1 We note that as we vary j3, anisotropics and 
masses, the values predicted by Eq. [40] agree to a few percent. 
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FIG. 4: The one-loop coefficient of the spatial plaquette, c sa , as a function of p for n p = 1, 2, 3. The minimum c aa (where the 
plaquette is near 1) in all three cases, as shown above, indicates p k, 0.14. 




Ensemble A: n p =2, p=0.22, P=2.0, Y g =3.5, m =0, y f =3.89 
o— o Ensemble B: n =1, p=0.14, [S=2.Q, y =3.5, m =0, y ( =3.5 
<^> Ensemble C: n =2, p=0.22, P=1.0, y =3.0, m =0, y t =3.33 



(1.2 
P 



FIG. 5: This graph shows how the spatial stout-smeared plaquette varies as a function of valence stout parameter, p, on three 
different ensembles with different choices of fermion parameters, (3 and stout-smearing factors, as indicated in the legend box. 
The smeared threshold is only observed for pthresh < 0.2 in j3 = 2.0 case; for smaller (3, the maximum /5 t hrcsn ~ 0.25. 



B. Gauge Anisotropy 



We determine the gauge anisotropy from the static-quark potential using Klassen's "Wilson-loop ratio" ap- 
proach [48|. In this method, we use Wilson loops involving the temporal direction, W s t, and those in the spatial 
directions, W ss , with Schrodingcr-functional boundary conditions applied in the z direction. We measure the ratios 

R ( T „\ - W ss( X >y) as y m , p -a 3 V 3 (ya 3 ) 

W ss {x + l,y) 



R st (x,t) 



Wst(?l tj asyni^ ^-a 3 V 3 {ta t ) 



W st (x + l,t) 



(42) 



on lattice volumes of 12 3 x 32. Naturally, one should impose R ss (x,y) — R s t(x,t) to get the renormalized £ g . An 
advantage to this method is that finite-volume artifacts tend to cancel in the ratios as demonstrated in the case of 
the (quenched) Wilson gauge action [48| and the Iwasaki gauge action [l?} with dynamical fermions. We determine 
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Self-consistent spatial tadpole: N f =3, m =0 
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FIG. 6: The self-consistent tadpole factors as a functions of j3 at fixed mo, 7 9 = 7/ ■ = £ — 3.5 and Ut = Ut = 1. The unsmcarcd 
spatial tadpole factors iti un ^ are nearly identical for both n p = 1 and n p = 2. 




2 2 
FIG. 7: Pade approximation for tadpole factors u s and u a for n p = 2. 



^ g by minimizing [l7j 



tic s _ {R ss (x,y) - R s t(x,£ g y)r 
" ^ (Ai? s )2 + (Ai? t ) 2 ' 



j'.y 



where Ai? s and Ai? t are the statistical errors of R ss and i? st . Wc interpolate R s t(x, t) by a cubic spline in terms of t. 
To avoid short-range lattice artifacts, x and y should not be too small. It has been observed in the Nf = 2 case [TtJ 
that including data with y = 1 introduces significant artifacts due to excited-state contamination. 

While this method was originally proposed for use in gauge actions with periodic boundary conditions, we note that 
it can also be applied to the case of a constant background field. The basic observation is that if we use Wilson loops 
orthogonal to the background field direction (the boundary condition direction), the flux from the background field 
that propagates should cancel in the ratios R ss or R st . To confirm this, we measured £ g with two different three-flavor 
ensembles with the same action parameters (3 = 1.5, j g = 4.4, mo — —0.0570 and 7/ = 3.3 but under two boundary 
conditions: Schrodinger-functional and periodic boundary conditions in the z direction. Figure [S] shows the results 
from both measurements. The gauge anisotropy from periodic boundaries, £g BC , as a function of mm(xy) is plotted. 
Compared with two-flavor simulation (l7| . the three-flavor results show larger dependence on miia(xy) . The band in 
Figure[5]is the resulting gauge anisotropy from Schrodinger-functional boundaries, £^ F , (obtained from min(xy) = 4); 
it is consistent with £j BC - We gain advantages by calculating £| F instead of £j BC for reasons described in Sec. Ill Cl 
Therefore, for the rest of this work, we will use the gauge anisotropy ^ F as £ g . 
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FIG. 8: Consistency check on £ 3 measured in Schrodinger-functional boundary conditions (band) and periodic boundary 
conditions (points) as functions of the product of x and y. The fermion action parameters for this run are 7 3 = 4.4, mo = —0.0570 
and 7/ = 3.3. 



C. Fermion Anisotropy 

We determine the fermion anisotropy £/ through the conventional relativistic meson dispersion relation 

E(j!) 2 = m 2 + £ (44) 

where the energy E and the mass m are in units of a t , and p = 1ixnjL Sl with L s the spatial lattice size, is in units 
of a s . Gauge configurations are generated using periodic boundary conditions in space and antiperiodic in time on 
lattice volumes of 12 3 x 96. We use six equally separated time sources, two different source quark smcarings and 
an unsmeared (local) sink to produce the correlators. These are then averaged over the six time sources. The two 
resulting hadron correlators for each momentum, averaged over equivalent rotations, were used in a constrained fit to 
two amplitudes and one mass. These results are also cross-checked against a fit including two masses (a ground and 
excited state). 

We calculate the energy E(p) at the spatial momenta p = 2irn/L s for n = (0,0,0), (1,0,0), (1,1,0), and (2,0,0). 
The resulting jackknife energies are used in a linear fit (in E(p) 2 ) to Eq. 03] to extract £/. We find no significant 
deviation from linearity. An example is in Figure [9] which shows the pseudoscalar and vector meson energies with 
action parameters /3 = 1.5, 7 9 ~ 4.4, m = —0.0540 and 7^ = 3.3 at multiple momentum projections. (The effective 
energy is shown as ln(C(i + 1)/C(t)).) We measure £/ = 3.44(8) from the pseudoscalar meson and £/ = 3.44(10) from 
the vector meson, which are consistent. 

In our earlier work, we determined the fermion anisotropy through the PCAC mass ratio from Schrodinger-functional 
boundary conditions in the z and t directions. However, unlike the case of the Wilson loop, we found it difficult to 
exclude excited-state contamination with the limited size of the spacial direction; thus, such a method suffers from 
a larger systematic error due to the extraction of the PCAC mass. We also investigated using the Schrodinger- 
functional scheme with Dirichlct (zero) gauge boundary conditions and a point-like (but smeared) fermion boundary 
condition 49] in time allowing us to project our correlators onto non-zero momenta, and thus allowing an extraction 
of the fermion anisotropy from the meson dispersion relation. This method, however, did not give as reliable results as 
in the antiperiodic case. Therefore, we choose to use the meson dispersion to determine £/ on (anti)pcriodic boundary 
gauge ensembles. 



D. Tuning Result 

With P fixed at 1.5, we only need to tune the parameters 7 9 , 7/ and mo. Our strategy is to choose some suitable 
estimates of the anisotropics, and make a coarse scan in mo to find the region of the critical point as determined 
from the PCAC mass M t , where, as described in Sec. Ill CI a background field in the t direction was used. Due to 
the nice properties of the algorithm using the background field, fairly short runs over roughly 1000 trajectories with 
measurements every fifth trajectory are sufficient for the purpose of determining M t . Given a mass regime of interest, 
the gauge anisotropy £ 9 is determined using a background in the z direction as described in Scc. lIVBl Again, roughly 
1000 trajectories are sufficient. The measurements of M t and £ s used 12 3 x 32 volumes, and provide a reasonable range 
of the bare parameters mo and 7 S for further determinations of the fermion anisotropy £/. These latter measurements 
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FIG. 9: Effective mass piots for various momentum projections of the pseudoscalar (left) and vector (right) mesons and 
corresponding fits. The fermion action parameters for this run are mo = —0.0540, Lt = 96, 7 9 = 4.4 and 7/ = 3.3 with 
(anti)periodic boundary conditions. 



used the method described in Sec. HVCl on 12 3 x 96 volumes with antiperiodic boundary conditions in time. Longer 
runs of roughly 2000 trajectories are used with measurements every fifth trajectory and binned twice, so they are 
effectively measured every tenth trajectory. 

Of course, each of these simulations is independent of the others. The ensemble parameters, and the results of the 
determinations of M t , £ g , and £/ are summarized in Table QJ Table ITT1 and Table \TU\ respectively. 

For each measurement, we fit the data according to the linear ansatz in Eq. 1161 The fit parameters are 

ao = 0.5(8), 01 = 0.41(14), a 2 = 0.39(8), o 3 = 1.3(7) 
b = -4.8(12), 6i = 0.6(2), b 2 = 1.74(17), b 3 = 4(3) 

c = 0.395(18), ci = -0.082(3), c 2 = 0.0194(17), c 3 = 1.282(16). (45) 

The values of % 2 per degree of freedom are 6.0/16, 13.2/8, 34.7/17 respectively. Figures [TU1 [TT1 and [T2l demonstrate 
how the fits work using two-dimensional slices of the parameter and measurement spaces. Both £ s and £/ show little 
sensitivity to the parameter too, and thus the slope is poorly determined in our fits. However, we do not expect this 
would affect our final result much since the changes in mo are O(10~ 2 ) and the £ g and £/ are O(l). The interactions 
of £ g and £/ also have large uncertainties in the fits; this is also expected since the running coefficients, 7/ and 
7 g , are around 3-4 and there is a long extrapolation to 0. In terms of fits to M t , too has the dominant effect, a 
few times larger than "f g while jf has one magnitude smaller contribution. We found that £ g and have positive 
linear dependence on the bare parameters, while for the case of M t versus 7 g , we found that M t increases while j g 
decreases. The 7 g parameter is the dominant factor in £ = a s /a t . Increasing it, a t will increase as well; this leads to 
the only non-dimensionless measurement M t becoming smaller in units of a t . Similar effects are observed when we 
parameterize the pion mass squared or rho-meson mass (which are measured on the antiperiodic boundary condition 
ensemble); see Figure [TBI and [T4l 

We can gain a qualitative understanding of the origin of the opposite sign in the slopes of j g versus 7/ for the mass 
measurements by considering the classical dispersion relation for the action in Eq. [5] For a small background chromo- 
electric field Fqi, we find that the lattice dispersion relation energy at zero spatial momentum E = 2sinh(atmo)/at 
satisfies for our choice of clover coefficient ct from Eq. [5] 

E 2 {1 + a t m ) = m 2 - (j^ + i ^ + sinh(a t m )^ ^ a 0i F 0i . (46) 
Thus, the derivatives dE 2 /d^ g oc — F^ and dE 2 /djf cx +Foi have opposite sign. 
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We impose the renormalization condition that in the chiral limit {£ g , M t } = {3.5, 3.5, 0}. By solving Eq.fTTlwith 
the parameters in Eq. 1451 we found that 

{m cri 7 ;, 7 ;} = {-0.080(6), 4.38(8), 3.44(7)}. (47) 

By fixing the renormalization condition in the general case, {£ s , £/, M t } — {3.5, 3.5, m q }, we obtain the desired action 
parameters as a functions of the bare quark mass, as shown in Fig. 1161 As emphasized in the discussions above, 7* 
and 7^ have very small quark mass dependence from the chiral limit up to the heaviest too used in this work, —0.057 
(which corresponds to m q ~ 0.03). Of course, too is linearly proportional to the bare quark mass. 

Figure [TS] shows a subset of the mo values with their corresponding gauge and fermion anisotropies. The two 
simulation points that arc consistent with the desired points of {£ ff ,£/} = {3.5,3.5} are consistent with Eq. 1471 

Based upon the negligible mass dependence that is observed, we henceforth fix the anisotropy bare parameters to 

7* = 4.3, 7} = 3.4 (48) 
which corresponds to a parameter set used in Figure [151 The corresponding clover coefficients are 

c s = 1.589, ct = 0.903. 

Figure [17] shows a PCAC measurement as a function of too with [3 = 1.5, 7ff = 4.4, and 7/ = 3.3; and linear 
interpolation to M t = gives us 

a t m cr = -0.0854(5) . (49) 

Finally, we will discuss how our initial tadpole-improved clover coefficients c s ,t differ from the nonperturbative 
coefficients in the Schrodinger-functional scheme (as described in Sec. Ill C|) . Table [T] lists the PCAC mass, AM t and 
its tree- level value (in units of a^ 1 ) measured with the i-direction Schrodinger-functional boundary condition. If 
we had nonperturbatively tuned c S) t and imposed the renormalization conditions from Eq. 1151 we should expect that 
AM t = AM t (0) . There are four sets of ensemble parameters that satisfy the conditions on 7* and 7^ . If we extrapolate 
AM t to TO cr , we find AM t = —0.00022(57), which is about 1.5 standard deviations away from the tree-level value: 
—0.00167. (See Figure [TBI) We conclude that the tadpole-corrected tree-level coefficients with stout-link smearing 
are close enough to the nonperturbative 0(a)-improved coefficients in the three-flavor dynamical simulation. 



V. CONCLUSION AND OUTLOOK 



This is the first calculation that combines Schrodinger functional with stout-link smearing on three- flavor anisotropic 
clover action. The Schrodinger-functional boundary conditions allow us to work on small volumes and small quark 
masses with improved signal, which improves our ability to probe the chiral limit. Stout-link smearing improves the 
chiral and scaling properties of our chirally broken fermion action. We take advantage of the much smaller cost by 
implementing both in our calculation. 

We studied a range of stout-link parameters that may be safely applied to multiple choices of (3 and different 
coefficients in the fermion sector in three-flavor simulation. Although our numerical simulations suggested that higher 
values of p could be applied, in the end, the stout-link parameters we chose conservatively to be p = 0.14 and n p = 2, 
as suggested by a one-loop pcrturbativc calculation. 

In a preliminary three-flavor study of various (3 values, we found through the static quark potential that (3 = 1.5 
gave us the desired spatial lattice spacing around 0.12 fm in the physical limit. The remaining coefficients in the 
gauge action and clover coefficients c s .t in the fermion action were set to their tree-level tadpole improved values, 
which were numerically determined and interpolated using a Pade approximation. 

We determined the gauge anisotropy using Wilson-loop ratios and Schrodinger-functional boundary conditions in 
the z direction; we determined the fermion anisotropy using antiperiodic boundary conditions in time. The quark 
mass was estimated from PCAC measurements using Schrodinger-functional boundary conditions in the t direction 
from which we found the critical mass. 

We then tuned the remaining parameters t; g j to achieve the desired renormalizcd anisotropy £ = 3.5 in the chiral 
limit. We found {m cr ,j*,jj} = {—0.080(6), 4.38(8), 3.44(7)} from a linear parameterization. The mass dependence 
of the gauge and fermion anisotropies is found to be negligible. Hence, we have for the final parameters 

7^ = 4.3, 7} = 3.4, c s = 1.589, c t = 0.903, 

with atm cr = —0.0854(5). Further, we showed in the Schrodinger-functional scheme that when using stout-link smear- 
ing and numerically determined tadpole factors, our fermion action automatically fulfills an (on-shell) O(a)-improvod 
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FIG. 10: Gauge anisotropy as a function of 7/, 7 9 at various fixed mo. The straight lines are the fit functions in Eq. 1 161 keeping 
the remaining two parameters fixed. Table ITT1 lists the details of the ensemble parameters. 
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FIG. 11: Fermion anisotropy as a function of 7/, 7 S at various fixed mo. Table ITTT1 lists the details of the ensemble parameters. 
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TABLE I: PCAC mass (in units of the inverse temporal lattice spacing) measured from the f-direction Schrodinger functional 
boundary condition ensemble with volume 12 3 x 32. 



m 


7s 


7/ 


M t 


AM, 


AA/ t (0) 


-0.0950 


4.3 


3.5 


-0.0122(9) 


0.0003(10) 


-0.001547 


-0.0950 


4.3 


3.4 


-0.0121(9) 


0.0003(9) 


-0.00167 


-0.0950 


4.3 


3.3 


-0.0141(8) 


0.0002(9) 


-0.001798 


-0.0734 


4.3 


3.5 


0.0160(9) 


-0.0007(9) 


-0.001547 


-0.0734 


4.3 


3.4 


0.0149(7) 


-0.0007(6) 


-0.00167 


-0.0734 


4.3 


3.3 


0.0139(7) 


0.0006(11) 


-0.001798 


-0.0618 


4.2 


3.5 


0.0431(11) 


0.0002(4) 


-0.001427 


-0.0618 


4.2 


3.4 


0.036(2) 


-0.0004(8) 


-0.001545 


-0.0618 


4.2 


3.3 


0.0339(11) 


0.0004(8) 


-0.001672 


-0.0618 


4.3 


3.5 


0.0321(9) 


0.0003(5) 


-0.001547 


-0.0618 


4.3 


3.4 


0.0303(6) 


-0.0001(5) 


-0.00167 


-0.0618 


4.3 


3.3 


0.0297(6) 


-0.0003(4) 


-0.001798 


-0.0618 


4.4 


3.4 


0.0218(5) 


-0.0004(5) 


-0.001791 


-0.0618 


4.4 


3.3 


0.0213(7) 


0.0001(6) 


-0.001923 


-0.0570 


4.3 


3.4 


0.0349(7) 


0.0002(5) 


-0.00167 


-0.0570 


4.3 


3.3 


0.0342(10) 


0.0004(7) 


-0.001798 


-0.0570 


4.3 


3.2 


0.0311(12) 


0.0005(8) 


-0.001935 


-0.0570 


4.3 


3.1 


0.0300(9) 


0.0025(8) 


-0.002078 


-0.0570 


4.4 


3.3 


0.0248(12) 


0.0002(8) 


-0.001923 


-0.0570 


4.4 


3.2 


0.0239(6) 


0.0005(7) 


-0.002065 


-0.0570 


4.4 


3.1 


0.0220(7) 


0.0012(11) 


-0.002212 



TABLE II: Renormalized gauge anisotropy measured from the z-direction Schrodinger functional boundary condition ensemble 
with volume 12 3 x 32. 



m 


7s 


7/ 




-0.0950 


4.3 


3.5 


3.48(5) 


-0.0950 


4.3 


3.4 


3.42(3) 


-0.0950 


4.3 


3.3 


3.40(3) 


-0.0743 


4.3 


3.4 


3.47(10) 


-0.0734 


4.3 


3.4 


3.46(4) 


-0.0618 


4.2 


3.5 


3.48(4) 


-0.0618 


4.2 


3.4 


3.41(3) 


-0.0618 


4.2 


3.3 


3.42(2) 


-0.0618 


4.3 


3.5 


3.50(4) 


-0.0618 


4.3 


3.4 


3.47(4) 


-0.0618 


4.3 


3.3 


3.43(4) 


-0.0618 


4.3 


3.2 


3.38(7) 


-0.0618 


4.4 


3.3 


3.47(5) 


-0.0570 


4.3 


3.4 


3.48(7) 


-0.0570 


4.3 


3.3 


3.43(9) 


-0.0570 


4.3 


3.2 


3.39(3) 


-0.0570 


4.3 


3.1 


3.36(4) 


-0.0570 


4.4 


3.3 


3.50(4) 


-0.0570 


4.4 


3.2 


3.54(5) 


-0.0570 


4.4 


3.1 


3.399(16) 
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FIG. 12: PCAC mass (in units of the temporal lattice spacing) as a function of jf, 7 9 at various fixed mo. Table [T] lists the 
details of the ensemble parameters. 



TABLE III: Fermion anisotropy, vector-meson mass and pseudoscalar mass squared (in units of the inverse temporal lattice 
spacing) from the periodic boundary condition ensemble with volume 12 s x 96. 



7g 7/ 



-0.0743 4.3 

-0.0618 4.2 

-0.0618 4.2 

-0.0618 4.2 

-0.0618 4.3 

-0.0618 4.4 

-0.0570 4.3 

-0.0570 4.3 

-0.0570 4.3 

-0.0570 4.4 

-0.0570 4.4 

-0.0570 4.4 



3.4 
3.5 
3.4 
3.3 
3.4 
3.3 
3.3 
3.2 
3.1 
3.3 
3.2 
3.1 



m n /m p 



3.43(4 
3.62(5 
3.38(4 
3.18(3 
3.47(6 
3.25(6 
3.23(4 
3.19(5 
2.99(4 
3.43(6 
3.22(8 
2.91(11 



0.1501(9) 
0.2830(9) 
0.2753(10) 
0.2604(10) 
0.2232(15) 
0.1639(17) 
0.2401(13) 
0.2290(16) 
0.2164(18) 
0.193(3) 
0.182(3) 
0.151(3) 



0.222(3 
0.348(2 
0.337(2 
0.319(5 
0.290(4 
0.217(5 
0.299(4 
0.292(4 
0.261(7 
0.255(6 
0.233(5 
0.194(6 



0.677(9) 
0.814(5) 
0.816(5) 
0.817(11) 
0.769(9) 
0.754(16) 
0.804(8) 
0.784(10) 
0.828(20) 
0.758(17) 
0.780(17) 
0.78(2) 



renormalization condition. In particular, the clover coefficients c Sj t are consistent with being nonperturbatively tuned. 
We will apply the same approach for future tuning as we move to finer lattices. 

With the determined coefficients from this work, we are currently generating Nf = 2 + 1-flavor ensembles with 
multiple masses. Further measurements are being made during the gauge generation process to precisely determine 
the lattice spacing in the chiral limit and some hadronic properties [44[ . 
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FIG. 13: Pion mass squared (in units of a t 2 ) as a function of 7/, 7 9 at various fixed mo. Table Hill lists the details of the 
ensemble parameters. 
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FIG. 14: Rho-meson mass (in units of a t 1 ) as a function of 7/, 7 9 at various fixed mo. Tablc[TTT]lists the details of the ensemble 
parameters. 
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FIG. 16: Tuned gauge and fermion action parameters, mg (top), 7g (middle) and 7^ (bottom), as functions of the bare quark 
mass, a t m q , corresponding to the solution of Eg. 1171 There is negligible mass dependence in the tuned anisotropy parameters. 
With a spatial lattice spacing a s = 0.12 fm, the maximum extent of the horizontal axis is about 175 MeV. 




FIG. 17: PCAC mass at j3 — 1.5, 7 9 = 4.3 and 7/ = 3.4 as a function of mo. A linear interpolation 
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FIG. 18: The measured AM t (down triangles) in units of a t as a function of mo at fixed 7 S = 4.3 and 7/ = 3.4, along with 
AM t (0) (dashed line) 
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